#ifndef I_LIB_IP_I_NMS_H
#define I_LIB_IP_I_NMS_H

#include <cassert>

#include "../core/i_basic.h"

namespace idl
{
    /*Fast approximated non-maximum suppression (anms) within a 2D circular window with diameter 5, the input "score" 2D (height x width) buffer is assumed to
    be sparse, i.e., not all entries are valid. This routine only considers entries whose scores are larger than "empty".
    For speed consideration, the 2-pixel width boundary pixels are ignored.
    Note this routine speeds up the nms process using some heuristic assumptions, thus the results are not equivalent to
    a standard implementation, although it works just fine in practice*/
    template <typename T>
    inline void i_approx_sparse_nms_circular_d5(T* scores, int width, int height, T empty)
    {
        int h1m2 = height - 2;
        int w1m2 = width - 2;
        int w2m1 = width + width - 1;
        int i, j, row, col, rw = width + width;
        T si;
        for (row = 2; row < h1m2; ++row)
        {
            for (col = 2; col < w1m2; ++col)
            {
                i = rw + col;
                si = scores[i];
                
                //ignore entries whose value == empty
                if (si <= empty)
                { 
                    continue;
                } 
                
                //1st:
                j = i + 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //2nd:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //3rd:
                j = i + w1m2;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //4th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //5th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //6th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //7th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //8th:
                j = i + w2m1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //9th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //10th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
            }
            rw += width;
        }
    }

  /*Fast approximated fast non-maximum suppression within a 2D circular window with diameter 5,
    the input "coords" is a nx2 array stores the inhomogeneous feature candidate coordinates (x, y, ...),
    the input "score" 2D (height x width) buffer is assumed to be sparse, i.e., not all entries are valid.
    This routine only considers entries whose scores are larger than "empty".
    Note: for speed consideration, the 2-pixel width boundary pixels are ignored. This requires that the input "coords"
    buffer should not contain any features that are out of the 2-pixel boundary, otherwise can casue memory segmentation fault problem.
    Note this routine speeds up the nms process using some heuristic assumptions, thus the results are not equivalent to
    a standard implementation, although it works just fine in practice*/
    template <typename T, typename S>
    inline void i_approx_sparse_nms_circular_d5(const S* coords, 
        int n, 
        T* scores, 
        S width, 
        S height, 
        T empty)
    {
        S h1m2 = height - 2;
        S w1m2 = width - 2;
        S w2m1 = width + width - 1;
        int n2 = n + n;
        int c;
        S row, col, i, j;
        T si;
        for (c = 0; c < n2; c += 2)
        {
            col = coords[c];
            row = coords[c + 1];

            if (col < 2 || col >= w1m2 || row < 2 || row >= h1m2)
            {
                continue;
            }

            i = row * width + col;
            si = scores[i];
            
            if (si <= empty)
            { 
                continue;
            } 
            //ignore entries whose value == empty

            //1st:
            j = i + 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //2nd:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //3rd:
            j = i + w1m2;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //4th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //5th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //6th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //7th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //8th:
            j = i + w2m1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //9th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //10th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
        }
    }

  /*Fast approximated fast non-maximum suppression within a 2D circular window with diameter 7, the input "score" 2D (height x width) buffer is assumed to
    be sparse, i.e., not all entries are valid. This routine only considers entries whose scores are larger than "empty".
    For speed consideration, the 3-pixel width boundary pixels are ignored.
    Note this routine speeds up the nms process using some heuristic assumptions, thus the results are not equivalent to
    a standard implementation, although it works just fine in practice*/
    template <typename T>
    inline void i_approx_sparse_nms_circular_d7(T* scores, 
        int width, 
        int height, 
        T empty)
    {
        int h1m3 = height - 3;
        int w1m3 = width - 3;
        int w2m2 = width + width - 2;
        int w3m1 = w2m2 + width + 1;
        int i, j, row, col, rw = width + width + width;
        T si;
        for (row = 3; row < h1m3; ++row)
        {
            for (col = 3; col < w1m3; ++col)
            {
                i = rw + col;
                si = scores[i];

                if (si <= empty)
                { 
                    continue;
                } //ignore entries whose value == empty

                //1st:
                j = i + 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //2nd:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //3rd:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //4th:
                j = i + w1m3;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //5th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //6th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //7th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //8th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //9th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //10th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //11th:
                j = i + w2m2;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //12th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //13th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //14th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //15th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //16th:
                j = i + w3m1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //17th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //18th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
            }
            rw += width;
        }
    }

    /*Fast approximated fast non-maximum suppression within a 2D circular window with diameter 7, 
    the input "coords" is a nx2 array stores the inhomogeneous feature candidate coordinates (x, y, ...),
    the input "score" 2D (height x width) buffer is assumed to be sparse, i.e., not all entries are valid. 
    This routine only considers entries whose scores are larger than "empty".
    Note: for speed consideration, the 3-pixel width boundary pixels are ignored. This requires that the input "coords"
    buffer should not contain any features that are out of the 3-pixel boundary, otherwise can casue memory segmentation fault problem.
    Note this routine speeds up the nms process using some heuristic assumptions, thus the results are not equivalent to
    a standard implementation, although it works just fine in practice*/
    template <typename T, typename S>
    inline void i_approx_sparse_nms_circular_d7(const S* coords, 
        int n, 
        T* scores, 
        S width, 
        S height, 
        T empty)
    {
        S h1m3 = height - 3;
        S w1m3 = width - 3;
        S w2m2 = width + width - 2;
        S w3m1 = w2m2 + width + 1;
        int n2 = n + n;
        int c;
        S row, col, i, j;
        T si;
        for (c = 0; c < n2; c += 2)
        {
            col = coords[c];
            row = coords[c + 1];
            
            if (col < 3 || col >= w1m3 || row < 3 || row >= h1m3)
            {
                continue;
            }
            
            i = row * width + col;
            si = scores[i];
            if (si <= empty)
            { 
                continue;
            } //ignore entries whose value == empty
            
            //1st:
            j = i + 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //2nd:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //3rd:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //4th:
            j = i + w1m3;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //5th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //6th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //7th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //8th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //9th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //10th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //11th:
            j = i + w2m2;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //12th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //13th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //14th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //15th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //16th:
            j = i + w3m1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //17th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //18th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
        }
    }
 
    /*Fast approximated fast non-maximum suppression within a 2D circular window with diameter 9, the input "score" 2D (height x width) buffer is assumed to
    be sparse, i.e., not all entries are valid. This routine only considers entries whose scores are larger than "empty".
    For speed consideration, the 4-pixel width boundary pixels are ignored.
    Note this routine speeds up the nms process using some heuristic assumptions, thus the results are not equivalent to
    a standard implementation, although it works just fine in practice*/
    template <typename T>
    inline void i_approx_sparse_nms_circular_d9(
        T* scores, 
        int width, 
        int height, 
        T empty)
    {
        int h1m4 = height - 4;
        int w1m4 = width - 4;
        
        int w2 = width + width;
        int w2m3 = w2 - 3;
        int w3m2 = w2 + width - 2;
        int w4m1 = w2 + w2 - 1;
        int i, j, row, col, rw = w2 + w2;
        T si;
        for (row = 4; row < h1m4; ++row)
        {
            for (col = 4; col < w1m4; ++col)
            {
                i = rw + col;
                si = scores[i];

                if (si <= empty)
                { 
                    continue; 
                } //ignore entries whose value == empty

                //1st:
                j = i + 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //2nd:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //3rd:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //4th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //5th:
                j = i + w1m4;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //6th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //7th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //8th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //9th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //10th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //11th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //12th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //13th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //14th:
                j = i + w2m3;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //15th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //16th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //17th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //18th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //19th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //20th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //21th:
                j = i + w3m2;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //22th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //23th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //24th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //25th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //26th:
                j = i + w4m1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //27th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
                //28th:
                j += 1;
                if (scores[j] > si)
                {
                    scores[i] = empty;
                    continue;
                }
                else
                {
                    scores[j] = empty;
                }
            }
            rw += width;
        }
    }

    /*Fast approximated fast non-maximum suppression within a 2D circular window with diameter 9,
    the input "coords" is a nx2 array stores the inhomogeneous feature candidate coordinates (x, y, ...),
    the input "score" 2D (height x width) buffer is assumed to be sparse, i.e., not all entries are valid.
    This routine only considers entries whose scores are larger than "empty".
    Note: for speed consideration, the 4-pixel width boundary pixels are ignored. This requires that the input "coords"
    buffer should not contain any features that are out of the 4-pixel boundary, otherwise can casue memory segmentation fault problem.
    Note this routine speeds up the nms process using some heuristic assumptions, thus the results are not equivalent to
    a standard implementation, although it works just fine in practice*/
    template <typename T, typename S>
    inline void i_approx_sparse_nms_circular_d9(
        const S* coords, 
        int n, 
        T* scores, 
        S width, 
        S height, 
        T empty)
    {
        S h1m4 = height - 4;
        S w1m4 = width - 4;
        S w2 = width + width;
        S w2m3 = w2 - 3;
        S w3m2 = w2 + width - 2;
        S w4m1 = w2 + w2 - 1;
        int n2 = n + n;
        int c;
        S row, col, i, j;
        T si;
    
        for (c = 0; c < n2; c += 2)
        {
            col = coords[c];
            row = coords[c + 1];

            if (col < 4 || col >= w1m4 || row < 4 || row >= h1m4)
            {
                continue;
            }

            i = row * width + col;
            si = scores[i];

            //ignore entries whose value == empty
            if (si <= empty)
            { 
                continue;
            } 

            //1st:
            j = i + 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //2nd:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //3rd:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //4th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //5th:
            j = i + w1m4;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //6th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //7th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //8th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //9th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //10th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //11th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //12th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //13th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //14th:
            j = i + w2m3;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //15th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //16th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //17th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //18th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //19th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //20th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //21th:
            j = i + w3m2;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //22th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //23th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //24th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //25th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //26th:
            j = i + w4m1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //27th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
            //28th:
            j += 1;
            if (scores[j] > si)
            {
                scores[i] = empty;
                continue;
            }
            else
            {
                scores[j] = empty;
            }
        }
    }


    /*Extract 5x5 local maxima in the region of s starting at [top,left] and ending at (bottom,right), 
    top and left should be at least 2, bottom should be at most height-2, right should be at most width-2 (the routine forces these boundary conditions).
    Put the results in "coords" (stores x0, y0, x1, y1, ...)  and return the number of local maximas.
    For safety, "coords" should be pre-allocated to the size of [2 x max_nr];
    Stop at max_nr (unsorted).
    Note: we assume that the abs operation has already been taken to elements in score buffer s, therefore s contains only positive values.*/
    template <typename T, typename S>
    inline int i_dense_nms_squared_5x5(const T *const* s, 
        S* coords,
        int width, int height,
        int left, int top, int right, int bottom, int max_nr,
        T abs_threshold)
    {
        int r, c, pos = 0;
        assert(left < right && top < bottom);
        /*actual boundary values*/
        int aleft = i_max(2, left);
        int atop  = i_max(2, top);
        int aright = i_max(width - 2, right);
        int abottom = i_max(height - 2, bottom);
        T v;
        for (r = atop; r < abottom; r++)
        {
            for (c = aleft; c < aright; c++)
            {
                v = s[r][c];
                /*25 comparisions*/
                if (v >= abs_threshold &&
                    v > s[r - 2][c - 2] && v > s[r - 2][c - 1] && v > s[r - 2][c] &&
                    v > s[r - 2][c + 1] && v > s[r - 2][c + 2] &&
                    v > s[r - 1][c - 2] && v > s[r - 1][c - 1] && v > s[r - 1][c] &&
                    v > s[r - 1][c + 1] && v > s[r - 1][c + 2] &&
                    v > s[r][c - 2] && v > s[r][c - 1] && v > s[r][c + 1] && v > s[r][c + 2] &&
                    v > s[r + 1][c - 2] && v > s[r + 1][c - 1] && v > s[r + 1][c] &&
                    v > s[r + 1][c + 1] && v > s[r + 1][c + 2] &&
                    v > s[r + 2][c - 2] && v > s[r + 2][c - 1] && v > s[r + 2][c] &&
                    v > s[r + 2][c + 1] && v > s[r + 2][c + 2])
                {
                    if (max_nr > (pos >> 1))/*divided by 2*/
                    {
                        /*feature*/
                        coords[pos++] = (S)c;
                        coords[pos++] = (S)r;
                    }
                }
            }
        }
        return (pos >> 1);/*divided by 2*/
    }

    /*Extract 5x5 local maxima in the region of s starting at [top,left] and ending at (bottom,right),
    top and left should be at least 2, bottom should be at most height-2, right should be at most width-2 (the routine forces these boundary conditions).
    Put the results (x, y, s) in xb, yb, and sb, and return the number of local maximas.
    For safety, "xb, yb, sb" should be pre-allocated to the size of "max_nr";
    Stop at max_nr (unsorted).
    Note: we assume that the abs operation has already been taken to elements in score buffer s, therefore s contains only positive values.*/
    template <typename T, typename S>
    inline int i_dense_nms_squared_5x5(const T *const* s,
        /*buffers that store coordinates x, y, and corner response value s, respectively*/
        S* xb, S* yb, T* sb,
        int width, int height,
        int left, int top, int right, int bottom, int max_nr,
        T abs_threshold)
    {
        int r, c, count = 0;
        //assert(left < right && top < bottom);
        /*actual boundary values*/
        int aleft = i_max(2, left);
        int atop = i_max(2, top);
        int aright = i_min(width - 2, right);
        int abottom = i_min(height - 2, bottom);
        T v;
        for (r = atop; r < abottom; r++)
        {
            for (c = aleft; c < aright; c++)
            {
                v = s[r][c];
                /*25 comparisions*/
                if (v >= abs_threshold &&
                    v > s[r - 2][c - 2] && v > s[r - 2][c - 1] && v > s[r - 2][c] &&
                    v > s[r - 2][c + 1] && v > s[r - 2][c + 2] &&
                    v > s[r - 1][c - 2] && v > s[r - 1][c - 1] && v > s[r - 1][c] && 
                    v > s[r - 1][c + 1] && v > s[r - 1][c + 2] &&
                    v > s[r][c - 2] && v > s[r][c - 1] && v > s[r][c + 1] && v > s[r][c + 2] &&
                    v > s[r + 1][c - 2] && v > s[r + 1][c - 1] && v > s[r + 1][c] &&
                    v > s[r + 1][c + 1] && v > s[r + 1][c + 2] &&
                    v > s[r + 2][c - 2] && v > s[r + 2][c - 1] && v > s[r + 2][c] &&
                    v > s[r + 2][c + 1] && v > s[r + 2][c + 2])
                {
                    if (max_nr > count)
                    {
                        /*feature*/
                        xb[count] = (S)c;
                        yb[count] = (S)r;
                        sb[count] = v;
                        ++count;
                    }
                }
            }
        }
        return (count);
    }


    /*Extract 5x5 local maxima in the region of s starting at [top,left] and ending at (bottom,right),
    top and left should be at least 2, bottom should be at most height-2, right should be at most width-2 (the routine forces these boundary conditions).
    Put the results (x, y, s) in xb, yb, and sb, and return the number of local maximas.
    For safety, "xb, yb, sb" should be pre-allocated to the size of "max_nr";
    Stop at max_nr (unsorted).
    Note: we assume that the abs operation has already been taken to elements in score buffer s, therefore s contains only positive values.*/
    inline int i_dense_nms_squared_5x5_f(const float *const* s,
        int* xb, 
        int* yb,
        /*buffers that store coordinates x, y, and corner response value s, respectively*/
        float* sb,
        int width, int height,
        int left, int top, int right, int bottom, int max_nr,
        float abs_threshold)
    {
        int r, c, count = 0;
        assert(left < right && top < bottom);
        /*actual boundary values*/
        int aleft = i_max((int)2, left);
        int atop = i_max((int)2, top);
        int aright = i_min(width - 2, right);
        int abottom = i_min(height - 2, bottom);
        float v;
        for (r = atop; r < abottom; r++)
        {
            for (c = aleft; c < aright; c++)
            {
                v = s[r][c];
                /*25 comparisions*/
                if (v >= abs_threshold &&
                    v>s[r - 2][c - 2] && v>s[r - 2][c - 1] && v > s[r - 2][c] &&
                    v > s[r - 2][c + 1] && v > s[r - 2][c + 2] &&
                    v > s[r - 1][c - 2] && v > s[r - 1][c - 1] && v > s[r - 1][c] &&
                    v > s[r - 1][c + 1] && v > s[r - 1][c + 2] &&
                    v > s[r][c - 2] && v > s[r][c - 1] && v > s[r][c + 1] && v > s[r][c + 2] &&
                    v > s[r + 1][c - 2] && v > s[r + 1][c - 1] &&
                    v > s[r + 1][c] && v > s[r + 1][c + 1] && v > s[r + 1][c + 2] &&
                    v > s[r + 2][c - 2] && v > s[r + 2][c - 1] &&
                    v > s[r + 2][c] && v > s[r + 2][c + 1] && v > s[r + 2][c + 2])
                {
                    if (max_nr > count)
                    {
                        /*feature*/
                        xb[count] = c;
                        yb[count] = r;
                        sb[count] = v;
                        ++count;
                    }
                }
            }
        }
        return (count);
    }

    /*Extract 5x5 local maxima in the region of s starting at [top,left] and ending at (bottom,right),
    top and left should be at least 2, bottom should be at most height-2, right should be at most width-2 (the routine forces these boundary conditions).
    Put the maxima values in sb and return the number of local maximas.
    "sb" should be pre-allocated to the size of "s".
    Note: we assume that the abs operation has already been taken to elements in score buffer s, therefore s contains only positive values.*/
    inline int i_dense_nms_squared_5x5_f(const float *const* s,
        /*buffer with non-maxima value in s are set to 0.0 in sb, maxima value set to the original response value in s*/
        float **sb, 
        int width, int height,
        int left, int top, 
        int right, int bottom,
        float abs_threshold)
    {
        int r, c, rm1, rm2, rp1, rp2, cm1, cm2, cp1, cp2;
        int count = 0;
        assert(left < right && top < bottom);
        
        memset((void*)sb[0], 0, width*height*sizeof(float));
        
        /*actual boundary values*/
        int aleft = i_max((int)2, left);
        int atop = i_max((int)2, top);
        int aright = i_min(width - 2, right);
        int abottom = i_min(height - 2, bottom);
        float v;
        for (r = atop; r < abottom; r++)
        {
            rm1 = r - 1;
            rm2 = r - 2;
            rp1 = r + 1;
            rp2 = r + 2;
            for (c = aleft; c < aright; c++)
            {
                v = s[r][c];
        
                cm1 = c - 1;
                cm2 = c - 2;
                cp1 = c + 1;
                cp2 = c + 2;
                
                /*25 comparisions*/
                if (v >= abs_threshold &&
                    v > s[rm2][cm2] && v > s[rm2][cm1] && v > s[rm2][c] && 
                    v > s[rm2][cp1] && v > s[rm2][cp2] &&
                    v > s[rm1][cm2] && v > s[rm1][cm1] && v > s[rm1][c] && 
                    v > s[rm1][cp1] && v > s[rm1][cp2] &&
                    v > s[r][cm2]   && v > s[r][cm1]   && v > s[r][cp1] && 
                    v > s[r][cp2] &&
                    v > s[rp1][cm2] && v > s[rp1][cm1] && v > s[rp1][c] && 
                    v > s[rp1][cp1] && v > s[rp1][cp2] &&
                    v > s[rp2][cm2] && v > s[rp2][cm1] && v > s[rp2][c] && 
                    v > s[rp2][cp1] && v > s[rp2][cp2])
                {
                    sb[r][c] = v;
                    ++count;
                }
            }
        }
        return (count);
    }

    /*Extract 5x5 local maxima in the region of s starting at [top,left] and ending at (bottom,right),
    top and left should be at least 2, bottom should be at most height-2, right should be at most width-2 (the routine forces these boundary conditions).
    Put the maxima values in sb and return the number of local maximas.
    "sb" should be pre-allocated to the size of "s".
    Note: we assume that the abs operation has already been taken to elements in score buffer s, therefore s contains only positive values.*/
    /*The 5x5 mask is rotated by 45 degrees so the actual number of comparisions is 12 instead of 24*/
    inline int i_dense_nms_squared_5x5_rot_f(const float *const* s,
        /*buffer with non-maxima value in s are set to 0.0 in sb, maxima value set to the original response value in s*/
        float **sb,
        int width, int height,
        int left, int top,
        int right, int bottom,
        float abs_threshold)
    {
        int r, c, rm1, rm2, rp1, rp2, cm1, cm2, cp1, cp2;
        int count = 0;
        assert(left < right && top < bottom);
    
        memset((void*)sb[0], 0, width*height*sizeof(float));
    
        /*actual boundary values*/
        int aleft = i_max((int)2, left);
        int atop = i_max((int)2, top);
        int aright = i_min(width - 2, right);
        int abottom = i_min(height - 2, bottom);
        float v;
        for (r = atop; r < abottom; r++)
        {
            rm1 = r - 1;
            rm2 = r - 2;
            rp1 = r + 1;
            rp2 = r + 2;
            for (c = aleft; c < aright; c++)
            {
                v = s[r][c];
    
                cm1 = c - 1;
                cm2 = c - 2;
                cp1 = c + 1;
                cp2 = c + 2;
    
                /*13 comparisions*/
                if (v >= abs_threshold &&
                    v > s[rm2][c] &&
                    v > s[rm1][cm1] && v > s[rm1][c] && v > s[rm1][cp1] &&
                    v > s[r][cm2] && v > s[r][cm1] && v > s[r][cp1] && v > s[r][cp2] &&
                    v > s[rp1][cm1] && v > s[rp1][c] && v > s[rp1][cp1] &&
                    v > s[rp2][c])
                {
                    sb[r][c] = v;
                    ++count;
                }
            }
        }
        return (count);
    }

    /*Extract 3x3 local maxima in the region of s starting at [top,left] and ending at (bottom,right),
    top and left should be at least 1, bottom should be at most height-1, right should be at most width-1 (the routine forces these boundary conditions).
    Put the maxima values in sb and return the number of local maximas.
    "sb" should be pre-allocated to the size of "s".
    Note: we assume that the abs operation has already been taken to elements in score buffer s, therefore s contains only positive values.*/
    inline int i_dense_nms_squared_3x3_f(const float *const* s,
        /*buffer with non-maxima value in s are set to 0.0 in sb, maxima value set to the original response value in s*/
        float **sb,
        int width, int height,
        int left, int top,
        int right, int bottom,
        float abs_threshold)
    {
        int r, c, rm1, rp1, cm1, cp1;
        int count = 0;
        assert(left < right && top < bottom);
        memset((void*)sb[0], 0, width*height*sizeof(float));

        /*actual boundary values*/
        int aleft = i_max((int)1, left);
        int atop = i_max((int)1, top);
        int aright = i_min(width - 1, right);
        int abottom = i_min(height - 1, bottom);
        float v;
        for (r = atop; r < abottom; r++)
        {
            rm1 = r - 1;
            rp1 = r + 1;
            for (c = aleft; c < aright; c++)
            {
                v = s[r][c];
                cm1 = c - 1;
                cp1 = c + 1;

                /*9 comparisions*/
                if (v >= abs_threshold &&
                    v > s[rm1][cm1] && v > s[rm1][c] && v > s[rm1][cp1] &&
                    v > s[r][cm1] && v > s[r][cp1] &&
                    v > s[rp1][cm1] && v > s[rp1][c] && v > s[rp1][cp1])
                {
                    sb[r][c] = v;
                    ++count;
                }
            }
        }
        return (count);
    }
} /* namespace idl */

#endif